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SECTION  1 


INTRODUCTION 

Seismic  event  location  remains  as  one  of  the  most  important  discriminants  for 
separating  natural  tectonic  and  explosive  events.  For  example,  it  has  been 
estimated  that  approximately  80%  of  known  global  earthquakes  have  focal  depths 
greater  than  50  km  or  are  located  more  than  25  km  at  sea.  Since  underwater 
nuclear  explosions  can  be  confidently  identified  to  very  low  yields  using 
hydroacoustic  data,  it  follows  that  the  vast  majority  of  earthquakes  can  potentially 
be  identified  on  the  basis  of  location  alone.  This  is  an  important  observation,  since 
global  monitoring  of  underground  nuclear  testing  may  require  the  screening  of  tens 
of  thousands  of  seismic  events  per  year  and,  in  order  to  perform  this  function 
efficiently,  it  will  be  necessary  to  have  simple  and  robust  discriminants  available 
which  can  be  used  to  eliminate  the  vast  majority  of  these  events  from  more 
detailed  consideration.  However,  in  order  to  be  useful  for  such  discrimination 
purposes,  the  uncertainties  associated  with  seismic  locations  must  be  well-defined 
and  reliable;  and  this  has  proven  to  be  difficult  to  accomplish  to  the  required 
degree  of  accuracy. 

More  specifically,  high  confidence  estimation  of  focal  depths  remains  as  an 
outstanding  seismic  monitoring  problem.  If  accurate  focal  depths  could  be 
determined  for  the  majority  of  earthquakes  deeper  than  5,  or  even  10  km,  the  event 
screening  problem  would  be  much  easier.  However,  despite  intensive  efforts  by 
well-trained  and  experienced  analysts,  more  than  two-thirds  of  the  GSETT-3  REB 
events  were  assigned  artificially  constrained  depths.  Moreover,  the  cited 
accuracies  of  even  those  relatively  few  focal  depth  estimates  is  open  to  serious 
question,  particularly  in  the  light  of  the  evidence  that  the  corresponding  epicentral 
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uncertainties  appear  to  be  significantly  underestimated  by  the  current  hypocentral 
location  algorithm.  Thus,  for  example,  Fisk  et  al  (2000)  found  a  number  of  cases  in 
the  GSETT-3  database  where  the  cited  uncertainties  in  the  focal  depths  of 
earthquakes  in  particular  areas  had  been  consistently  underestimated  by  the  PIDC 
location  program.  Such  discrepancies  are  frequently  compensated  for  in  a  crude 
fashion  by  artificially  inflating  the  input  variance  values  until  the  associated 
confidence  bounds  encompass  the  required  percentage  of  verified  values.  This  was 
the  interim  approach  taken  by  Fisk  et  al  (2000)  in  their  formulation  of  the  focal 
depth  discriminant  for  the  experimental  event  screening  system  at  the  IDC. 
However,  this  is  an  ad  hoc  procedure  which  can  potentially  lead  to  serious  focal 
depth  estimation  errors  under  previously  untested  conditions. 

This  lack  of  a  reliable  measure  of  focal  depth  uncertainty  has  a  number  of 
important  ramifications  with  respect  to  event  screening.  That  is,  not  only  does  it 
reduce  the  effectiveness  of  the  focal  depth  discriminant  itself,  but  it  also  has 
implications  for  some  of  the  other  standard  discriminants.  This  follows  from  the 
fact  that  from  a  seismic  monitoring  perspective,  if  an  event  can’t  be  confidently 
shown  to  have  significant  focal  depth,  it  must  be  assumed  to  be  shallow  enough  to 
be  considered  as  a  possible  explosion.  For  this  reason,  the  default  focal  depth 
typically  employed  in  seismic  monitoring  of  explosions  is  zero,  as  opposed  to  the 
33  km  value  frequently  employed  by  NEIS.  Now  if  events  which  are  actually 
subcrustal  are  assigned  shallow  focal  depths,  this  can  have  a  pronounced  effect  on 
the  utility  of  other  discriminants,  such  as  Mj/mb.  For  example,  it  has  been  shown 
that  deep,  subcrustal  earthquakes  generally  show  explosion-like  Ms-mb  values  (i.e., 
Ms-iUb  <  -1.0),  in  accord  with  the  theoretically  expected  decrease  in  Ms  with 
increasing  focal  depth.  It  follows  that,  if  events  which  are  actually  subcrustal  have 
to  be  analyzed  under  the  assumption  that  they  are  shallow,  the  effectiveness  of  the 
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Ms/mb  discriminant  will  be  greatly  reduced.  Thus,  any  improvements  that  can  be 
made  in  focal  depth  determination  capability  will  have  significant  and  multiple 
implications  with  respect  to  operational  monitoring  capability. 

This  report  presents  a  summary  of  a  variety  of  investigations  conducted  under  this 
contract  that  have  been  directed  toward  the  development  of  improved  procedures 
for  estimating  focal  depths  of  seismic  sources.  The  utility  of  the  secondary  depth 
phases  pP  and  sP  is  addressed  in  Section  2,  where  the  various  factors  affecting 
depth  phase  identification  are  discussed  and  illustrated;  and  a  network 
beamforming  procedure  for  use  in  the  automatic  identification  of  candidate  depth 
phases  is  described  and  applied  to  seismic  data  recorded  from  large  samples  of 
Hindu  Kush,  Honshu,  Hokkaido  and  Lop  Nor  earthquakes.  In  Section  3  we  discuss 
the  use  of  the  F  detector  as  a  tool  to  assist  in  determining  the  validity  of  candidate 
depth  phases.  This  is  followed  in  Section  4  with  an  assessment  of  the  use  of 
calibrated  S-P  differential  arrival  time  information  to  provide  independent 
estimates  of  event  origin  time  for  use  as  constraints  in  the  hypocenter  location 
process.  The  report  concludes  with  Section  5,  which  contains  a  short  summary  and 
statement  of  conclusions  regarding  improvements  in  focal  depth  determination  for 
use  in  the  seismic  monitoring  of  underground  nuclear  explosions. 
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SECTION  2 


IMPROVED  DETECTION  AND  IDENTIFICATION 
OF  SECONDARY  DEPTH  PHASES 


2.1  BACKGROUND. 

Reliable  identification  of  the  depth  phases  pP  and  sP  provides  one  of  the  strongest 
pieces  of  evidence  that  a  seismic  event  is  a  naturally  occurring  earthquake  when 
the  time  delay  between  pP  and  direct  P  exceeds  about  three  seconds,  implying  a 
focal  depth  of  greater  than  10  km.  However,  confident  identification  of  such 
phases  on  actual  recordings  is  not  always  a  simple  matter.  That  is,  it  is  first 
necessary  to  unambiguously  identify  a  phase  onset  in  the  P  coda  background;  and 
it  is  then  required  to  confirm  that  this  phase  is  actually  a  depth  phase  and  not  some 
other  P  phase,  such  as  PcP,  which  might  have  an  expected  arrival  time  within  the 
time  window  being  analyzed.  Thus,  with  regard  to  secondary  phase  detection, 
analysts  are  generally  required  to  verify  that  the  phase  amplitude  is  several  times 
the  P  coda  background  before  it  can  be  accepted  for  consideration  as  a  possible 
depth  phase.  For  positive  identification  as  a  depth  phase,  it  is  generally  further 
required  that  the  subject  arrival  be  detected  at  several  stations  at  different 
epicentral  distances  and  that  the  time  interval  between  its  onset  time  and  the 
associated  P  onset  time  be  an  increasing  function  of  epicentral  distance.  This 
distance  dependence  is  commonly  referred  to  as  “step-out,”  and  the  lASPEI  91 
travel  time  tables  predict  step-outs  between  30  and  90  degrees  of  approximately 
0.4  seconds  for  a  focal  depth  of  15  km,  4.4  seconds  for  a  focal  depth  of  100  km  and 
13.7  seconds  for  a  focal  depth  of  300  km.  Thus,  even  with  high  quality  data  for 
which  it  is  possible  to  determine  pP-P  time  of  arrival  differences  with  an  accuracy 
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of  about  1  second,  it  is  generally  not  possible  to  apply  the  step-out  identification 
criterion  to  pP  for  events  with  focal  depth  less  than  about  50  km.  Unfortunately, 
it  was  discovered  during  the  course  of  the  GSETT-3  that  it  is  not  an  easy  matter  to 
specify  analyst  review  procedures  which  will  assure  such  high  confidence 
identification  of  depth  phases  on  a  routine  basis.  For  example.  Figure  2-1  shows 
ARCES  and  NORES  data  recorded  from  the  Hindu  Kush  earthquake  of 
1997/02/22,  with  annotated  analyst  pP  picks  which  were  used  to  estimate  a  very 
precise  focal  depth  for  this  event.  It  can  be  seen  that,  while  there  is  evidence  of 
energy  arriving  near  the  picked  pP  times  at  both  these  stations,  the  arrivals  do  not 
appear  to  be  distinctively  different  from  the  preceding  P  coda;  and  their 
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Figure  2-1.  Center  element  seismograms  (sz)  and  infinite  velocity  beams 
(beam)  corresponding  to  the  ARCES  and  NORES  recordings  of  the 
1997/02/22  Hindu  Kush  earthquake  showing  the  analyst  picked  arrival  times 
for  P  and  pP. 
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identification  would  have  to  be  considered  to  be  highly  questionable  in  the 
absence  of  further  corroborative  evidence,  such  as  moveout.  However,  the 
moveout  criterion  can  not  be  applied  in  this  case  due  to  the  proximity  of  the  two 
stations;  and,  consequently,  these  analyst  picks  would  have  to  be  rejected  as 
unconfirmed  for  event  screening  purposes.  A  similar  example  is  shown  in 
Figure  2-2  for  the  Hindu  Kush  earthquake  of  1997/09/11.  In  this  case,  there  are 
three  reported  pP  observations  over  the  epicentral  distance  range  extending  from 
18  to  73  degrees,  which  do  show  some  indication  of  moveout.  However,  once 
again,  identification  of  the  picked  phases  on  the  recorded  waveforms  appears  to  be 
problematic  in  that  there  are  no  obvious  distinct  arrivals  at  the  analyst  picked  pP 
times  at  any  of  the  three  stations.  Experimentation  with  different  frequency  filters 


Figure  2-2.  Single  component  seismograms  corresponding  to  the  ZAL,  ARCES 
and  ILAR  recordings  of  the  1997/09/11  Hindu  Kush  earthquake  showing  the 
analyst  picked  arrival  times  for  P  and  pP. 
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and  array  beams  produced  no  better  results,  leading  to  the  conclusion  that  these  are 
also  false  detections.  There  are  also  numerous  examples  in  the  GSETT-3  database 
of  the  converse  problem,  in  which  the  analyst  failed  to  identify  fairly  obvious 
depth  phases,  leading  to  inferior  estimates  of  focal  depth.  Such  examples  indicate 
that  improved  tools  are  needed  to  help  the  analyst  to  improve  on  the  detection  and 
identification  of  the  depth  phases  pP  and  sP. 

Leaving  aside  for  the  moment  the  various  types  of  analyst  errors  that  can  occur,  it 
is  important  to  address  the  more  general  question  of  why  depth  phase  identification 
continues  to  be  such  a  difficult  problem  despite  all  the  advances  in  instmmentation 
and  signal  processing  which  have  been  implemented  in  recent  years.  For  example, 
during  the  GSETT-3  experiment,  approximately  2/3  of  the  events  reported  in  the 
REB  had  depths  artificially  constrained  to  the  surface;  and,  of  the  remaining  1/3  of 
the  events  which  were  assigned  free  depths,  only  about  20%  were  determined  on 
the  basis  of  depth  phase  data.  That  is,  depth  phase  data  were  reported  for  only 
about  5%  of  the  total  number  of  REB  events.  There  are  a  number  of  factors  which 
contribute  to  these  low  levels  of  detectability.  An  obvious  one  is  that,  for  the  many 
shallow  events,  the  time  separation  between  pP  and  P  may  be  so  short  in  some 
cases  that  the  pP  arrival  is  obscured  by  the  P  coda,  particularly  for  earthquakes 
with  complex  source  functions.  However,  the  problem  persists  for  deeper  events 
for  a  variety  of  more  subtle  reasons.  We  will  illustrate  some  of  these  complicating 
effects  using  data  recorded  from  REB  events  in  the  Hindu  Kush  region  of  Asia. 
Because  of  the  tectonics  in  this  region,  numerous  events  in  the  depth  range  200  to 
250  km  occur  in  close  proximity  to  one  another,  providing  a  convenient  data  set  to 
illustrate  depth  phase  variability.  The  first  effect  to  be  considered  is  associated 
with  the  fact  that  the  direct  P  and  surface  reflected  pP  and  sP  phases  do  not  travel 
exactly  the  same  paths  from  the  source  to  the  receiver;  and,  consequently,  their 
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pulse  shapes  can  appear  to  be  quite  different.  More  specifically,  the  surface 
reflected  phases  propagate  up  to  the  surface  and  back  through  the  attenuative 
surface  layers  above  the  source,  which  results  in  surface  reflected  pulses  which  are 
typically  of  lower  dominant  frequency  than  the  associated  direct  P  pulses.  This 
effect  is  illustrated  in  Figure  2-3  using  data  recorded  at  three  nearby  North 
American  stations  from  the  Hindu  Kush  earthquake  of  1997/01/10.  Note  that  the 
dominant  frequency  content  of  pP  is  lower  than  that  of  P  on  the  broadband  signals 
(left)  and  that  the  pP/P  amplitude  ratios  show  significant  variability  among  these 
three  stations  located  in  the  same  narrow  azimuth  window.  However,  when  these 
same  signals  are  lowpass  filtered  below  1.5  Hz  (right),  it  can  be  seen  that  the  pP/P 
amplitude  ratios  become  much  more  similar,  confirming  the  fact  that  the  large 
variability  observed  on  the  broadband  traces  is  due  to  differences  in  dominant 
frequency  content  between  the  P  and  pP  signals  observed  at  these  stations.  Such 
differences  contribute  to  the  difficulty  in  correlating  depth  phase  arrivals  between 
different  stations.  A  second  factor  which  complicates  the  depth  phase  identification 
process  is  earthquake  focal  mechanism  which  can  also  produce  large  variations  in 
signal  characteristics  between  stations  for  a  given  event  and  between  events  at  a 
given  station.  Such  effects  are  illustrated  in  Figure  2-4  which  shows  bandpass 
filtered  (1-2  Hz)  recordings  for  station  MBC  in  Canada  from  a  series  of  nearby 
Hindu  Kush  earthquakes  with  depths  in  the  200  to  250  km  range.  Note  the  order  of 
magnitude  variations  in  the  pP/P  and  sP/P  amplitude  ratios  which  are  presumably 
associated  with  variations  in  focal  mechanism  between  these  closely  spaced  events. 
More  surprising  is  the  fact  that  such  variability  in  focal  mechanism  between  nearby 
events  can  be  large  enough  in  some  cases  to  produce  consistent  differences  in  pP/P 
amplitude  ratios  at  large  numbers  of  widely  separated  stations.  An  example  of  this 
type  of  variability  is  shown  in  Figure  2-5  using  recordings  at  the  same  stations 
from  two  nearby  Hindu  Kush  earthquakes  with  comparable  focal  depth  (i.e.  208 
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Figure  2-3.  Comparison  of  broadband  (left,  0.6  -  4.5  Hz)  and  lowpass  filtered  (right,  0.8 
1.5  Hz)  recordings  of  the  southern  Hindu  Kush  event  of  1997/01/10  (h  =  207  km)  at  the 
three  nearby  North  American  Stations  MBC,  ILAR  and  YKA. 


h=237km 


h=240km 


h=248km 


10  sec 


Figure  2-4.  Comparison  of  station  MBC  (A=70°)  bandpass  filtered  (1-2  Hz) 
recordings  of  selected  Hindu  Kush  earthcjuakes  with  REB  depths  in  the  range 
200  <  h  <  250  km. 
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Figure  2-5.  Comparison  of  recordings  at  conmion  teleseismic  stations  from  two 
nearby  Hindu  Kush  events  of  comparable  magnitude  and  depth. 


and  228  km).  It  can  be  seen  that  for  one  of  these  earthquakes  (1998/12/11,  right) 
clear  depth  phases  are  evident  at  virtually  every  station;  whereas  for  the  other 
(1999/02/09,  left)  no  depth  phases  are  evident,  despite  the  fact  that  these  stations 
sample  an  azimuth  window  extending  over  nearly  180  degrees.  Again,  this  type  of 
variability  makes  it  difficult  to  correlate  arrivals  between  events,  even  in  narrowly 
confined  source  regions. 

A  third  complicating  factor  is  the  occasional  occurrence  of  anomalous  arrivals 
which  are  not  predicted  by  the  global  travel  time  tables.  Such  “mystery  phases” 
can  have  significant  amplitudes  and  could  be  easily  misidentified  as  depth  phases 
in  some  cases.  One  interesting  example  of  such  a  mystery  phase  is  shown  in 


10 


Figure  2-6  where  recordings  at  station  ESDC  in  Spain  from  a  number  of  nearby 
Hindu  Kush  earthquakes  are  displayed  in  order  of  increasing  focal  depth  between 
190  and  243  km.  Note  the  strong  anomalous  P  arrival  (P?)  in  the  interval  between 


P  and  pP  which  is  only  observed  over  a  limited  range  of  focal  depth.  In  the 
absence  of  other  information,  it  would  be  difficult  for  the  analyst  to  resist 
identifying  this  arrival  as  a  candidate  depth  phase,  particularly  in  those  cases  where 
its  amplitude  is  as  large  or  larger  than  that  of  the  true  later  arriving  pP  phase. 
Further  investigation  has  revealed  that  this  arrival  is  only  observed  from 
earthquakes  in  a  very  confined  focal  region.  Other  anomalous  P  phases  have  been 
observed  from  events  in  a  variety  of  different  source  regions  (eg.  Kind  and  Seidl, 
1982),  which  indicates  that  great  care  must  be  taken  to  conclusively  validate 
candidate  depth  phase  readings  before  they  can  be  used  to  conclude  with  high 
confidence  that  an  event  is  too  deep  to  be  an  explosion. 


Figure  2-6.  Comparison  of  station  ESDC  (Sonseca,  Spain)  recordings  of  selected 
Hindu  Kush  earthquakes  with  REB  depths  in  the  range  200  <  h  <  250  km. 


11 


2.2  APPLICATION  OF  NETWORK  BEAMFORMING  TO  THE  AUTOMATIC 
IDENTIHCATION  OF  CANDIDATE  DEPTH  PHASES. 

In  an  attempt  to  overcome  the  various  difficulties  in  single  station  depth  phase 
identification  illustrated  in  the  preceding  section,  we  have  been  investigating 
(Murphy  et  al,  2000)  the  potential  utility  of  an  automated  network  stacking 
algorithm  which  employs  signal  analysis  procedures  similar  to  those  originally 
proposed  by  Israelsson  (1994)  and  more  recently  applied  to  a  sample  of  Canadian 
data  by  Woodgold  (1999).  The  general  problem  is  represented  in  simplified  form 
in  Figure  2-7  which  shows  synthetic  pP  and  sP  arrivals  plotted  at  their  predicted 
delay  times  with  respect  to  P  for  a  hypothetical  Hindu  Kush  earthquake  with  a 
depth  of  150  km.  For  simplicity,  the  arrivals  are  represented  here  by  unit 
amplitude  boxcar  functions  with  constant  widths  of  ±1  second  for  pP  and  ±1.5 
second  for  sP.  The  trace  at  the  bottom  of  this  figure  shows  the  direct  sum  of  these 
individual  synthetic  traces,  and  it  can  be  seen  that  the  result  is  significantly 
degraded  because  the  moveouts  of  the  pP  and  sP  pulses  with  respect  to  P  as  a 
function  of  epicentral  distance  have  not  been  taken  into  account.  Israelsson  (1994) 
proposed  to  solve  this  problem  by  mapping  the  individual  post-P  arrivals  observed 
at  each  station  from  functions  of  delay  time  with  respect  to  P  to  functions  of  source 
depth  using  the  pP-P  or  sP-P  delay  times  predicted  by  the  lASPEI  91  travel  time 
tables  for  the  individual  station  distances.  This  one  to  one  mapping  is  illustrated  in 
Figure  2-8  which  shows  the  predicted  pP  moveout  with  respect  to  P  as  a  function 
of  source  depth  for  two  different  epicentral  distances.  Using  such  predicted 
moveout  curves,  the  observed  delay  times  of  all  post-P  arrivals  with  respect  to  P 
can  be  translated  into  equivalent  focal  depths  under  the  hypothesis  that  they  are  pP 
or  sP  arrivals.  Those  arrivals  which  are  consistent  with  such  hypotheses  should 
then  add  coherently  and  show  up  as  pronounced  peaks  on  the  stacked  network 
depth  function.  Figure  2-9  shows  the  result  of  applying  this  transformation  to  the 
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Figure  2-7.  Synthetic  pP  and  sP  arrivals  plotted  at  their  predicted  delay  times  with 
respect  to  P  for  a  hypothe<^ical  Hindu  Kush  earthquake  with  a  depth  of  150  km. 
The  bottom  trace  shows  the  direct  sum  of  these  individual  synthetic  arrivals. 


Figure  2-8.  pP-P  moveout  times  as  a  function  of  focal  depth  at  fixed  epicentral 
distances  of  30  and  80  degrees. 
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synthetic  arrivals  of  Figure  2-7  under  the  pP  hypothesis.  It  can  be  seen  that  the 
network  stack  of  these  transformed  synthetic  traces  produces  a  much  sharper  and 
larger  amplitude  peak  than  that  obtained  from  the  uncorrected  direct  sum  shown  in 
Figure  2-7. 


Figure  2-9.  Network  detection  stack  of  pP  as  a  function  of  candidate  source  depth 
corresponding  to  the  previous  synthetic  example. 


As  was  indicated  by  the  examples  shown  in  the  previous  section,  one  of  the 
difficulties  associated  with  the  application  of  the  above  algorithm  directly  to  the 
recorded  short-period  seismic  data  is  that  the  depth  phase  waveforms  and  dominant 
frequencies  for  a  given  event  can  show  substantial  variability  across  a  global 
network  due  to  both  focal  mechanism  and  propagation  path  effects,  which  would 
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significantly  degrade  the  effective  gain  of  any  direct  network  stacking  procedure. 
Consequently,  Israelsson  (1994)  and  Woodgold  (1999)  formulated  special  purpose 
signal  processing  algorithms  and  station  selection  criteria  in  attempts  to 
compensate  for  this  lack  of  signal  coherence  between  stations.  In  Murphy  et  al 
(2000)  an  alternate  approach  was  adopted  which  is  based  on  the  automatic 
detection  data  produced  by  the  detection  processing  employed  at  both  the  U.S. 
NDC  and  the  IDC.  In  this  approach,  the  observed  post-P  detection  times  in  the 
detection  file  for  a  given  event  are  used,  together  with  their  estimated  uncertainties 
(currently  taken  to  be  ±1  second  for  pP  and  ±1.5  seconds  for  sP),  to  define  unit 
amplitude  boxcar  functions  centered  on  the  automatic  detection  time;  and  these 
boxcars  as  functions  of  time  are  then  mapped  into  equivalent  boxcars  as  a  function 
of  source  depth  using  predicted  lASPEI  91  moveout  curves  such  as  those  shown  in 
Figure  2-8.  It  is  these  simplified  functions  of  depth  which  are  then  stacked  over 
the  network  of  observing  stations  to  identify  candidate  depth  phases  which  are 
consistent  with  either  the  pP  or  sP  hypothesis.  An  advantage  of  this  approach  is 
that  not  only  is  it  fully  automatic,  but  it  only  includes  signals  which  have  triggered 
one  of  the  detectors,  which  should  help  to  minimize  the  inclusion  of  questionable 
arrivals  in  the  analysis.  Figure  2-10  shows  an  example  of  the  application  of  this 
algorithm  to  the  automatic  post-P  detection  times  for  a  Hindu  Kush  earthquake  on 
1998/02/14.  It  can  be  seen  that  the  network  detection  stack  in  this  case  shows  a 
pronounced  peak  at  a  depth  corresponding  to  the  published  REB  depth  for  this 
event,  indicating  the  presence  of  a  good  candidate  pP  phase  which  should  be 
further  reviewed  by  the  analyst. 

In  Murphy  et  al  (2000)  the  depth  phase  stacking  algorithm  described  above  was 
first  applied  to  the  automatic  detection  data  from  over  150  selected  REB  events 


15 


0  100  200  300  400  500  600  700 


Depth,  km 

Figure  2-10.  Network  detection  stack  of  pP  for  the  Hindu  Kush  earthquake  of 
1998/02/14,  nib  =  5.03.  The  dashed  vertical  line  coincides  with  the  REB  depth 
estimate  of  226  km. 

located  in  the  Hindu  Kush,  Lop  Nor,  Hokkaido  and  central  Honshu  regions 
identified  on  the  map  of  Figure  2-11.  Selected  examples  are  shown  in  Figures  2-12 
through  2-15.  It  can  be  seen  that  in  all  these  examples  there  are  prominent  peaks 
on  the  network  detection  stacks  which  are  either  at  or  close  to  the  corresponding 
REB  depth  estimates.  In  some  cases,  the  REB  depths  do  show  some  measureable 
offset  with  respect  to  the  maxima  of  the  corresponding  network  detection  stacks. 
This  reflects  the  fact  that  no  depth  phase  data  were  used  in  determining  the  REB 
hypocentral  solutions  for  these  events  and  suggests  that  more  accurate  depths 
might  have  been  determined  for  these  events  if  the  candidate  depth  phases 
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Figure  2-11.  Map  locations  of  the  events  used  in  the  evaluation  of  the  depth  phase 
stacking  algorithm. 
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Figure  2-12.  Network  detection  stacks  of  pP  for  the  Hindu  Kush  earthquakes  of 
1998/02/14  (mb  =  5.03),  1995/07/09  (mb  =  4.40)  and  1996/04/08  (mb  =  4.04). 
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Figure  2-13.  Network  detection  stacks  of  pP  for  the  Hokkaido  earthquakes  of 
1997/03/13  (nib  =  4.34),  1997/07/22  (nib  =  3.97)  and  1996/04/13  (nib  =  3.82). 
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Figure  2-14.  Network  detection  stacks  of  pP  for  the  Central  Honshu  earthquakes 
of  1996/11/23  (mb  =  4.32),  1996/05/08  (mb  =  3.86)  and  1996/10/23  (mb  =  3.79). 
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Figure  2-15.  Network  detection  stacks  of  pP  for  the  Lop  Nor  earthquakes  of 
2001/03/13  (mb  =  4.11)  and  1999/01/27  (mb  =  4.01). 


identified  by  the  network  stacking  algorithm  could  have  been  validated  and 
included  in  the  hypocentral  solutions.  The  results  shown  in  Figures  2-12  through 
2-15,  along  with  other  results  we  obtained  (Murphy  et  al,  2000),  document  the  fact 
that  the  proposed  algorithm  is  applicable  to  events  having  mb  values  as  low  as 
about  4.0  in  the  Hindu  Kush  region  and  about  3.75  in  both  the  Hokkaido  and 
central  Honshu  regions  of  Japan,  witli  corresponding  minimum  depths  of  about  75 
km  and  55  km,  respectively.  Comparable  results  were  obtained  for  over  70  %  of 
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the  events  analyzed,  which  suggests  this  fully  automatic  algorithm  could  provide  a 
powerful  tool  for  directing  an  analyst  to  candidate  depth  phases  which  could  be 
used  in  the  event  location  process. 

2.3  INTEGRATION  INTO  THE  ANALYST  REVIEW  STATION 
The  Research  and  Development  Support  System  (RDSS)  developed  by  the  Army 
Space  and  Missile  Defense  Center  (SMDC)  was  established  to  improve  nuclear 
explosion  monitoring  capability  by  supporting  the  R+D  community  with  a  broad 
range  of  activities  and  resources  (Woodward  and  North,  2002).  RDSS  activities 
include  providing  environments  for  testing  and  evaluating  promising  research 
results  at  a  range  of  scales.  In  order  to  make  our  depth  phase  stacking  procedure 
available  for  large-scale  testing  and  further  evaluation  the  algorithm  described 
above  has  been  formally  integrated  into  the  RDSS  Analyst  Review  Station  (ARS) 
software  which  is  employed  at  both  the  U.S.  NDC  and  the  IDC.  This  milestone 
now  gives  us  the  capability  to  evaluate  our  algorithm  using  the  full-processing 
environment  of  the  representative  nuclear  monitoring  systems  as  well  as  historical 
data  archives.  During  the  interactive  processing  of  an  event,  the  analyst  is  able  to 
initiate  an  automatic  process  from  the  ARS  that  causes  all  of  the  post-P  detection 
times  for  the  current  event  of  interest  to  be  recovered  from  the  automatic  detection 
files  and  passed  to  the  depth  phase  stacking  module,  where  they  are  converted  into 
their  corresponding  depth  traces,  stacked,  and  presented  back  to  the  analyst  in  an 
interactive  X-Window  display.  This  interface  to  the  depth-stacking  algorithm  in 
ARS  is  illustrated  in  Figure  2-16  and  Figure  2-17,  where  we  first  show  an  example 
of  the  pP  network  detection  stack  for  a  Lop  Nor  region  earthquake  on  1999/01/27, 
followed  by  a  summary  of  the  pP,  sP,  and  pP  +  sP  stacks  for  this  event.  Note  that 
the  sum  of  the  pP  and  sP  staeks  shows  a  pronounced  candidate  depth  phase  peak 
near  the  published,  shallow  REB  depth  of  about  18  km  for  this  earthquake.  This 
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example  is  evidence  that,  at  least  in  some  cases,  this  depth  stacking  procedure  can 
be  applied  to  crustal  events  with  depths  shallower  than  50  km.  At  this  point  in  the 
ARS  implementation  of  the  procedure,  by  clicking  on  this  or  any  other  candidate 
depth-phase  peak  with  the  mouse,  the  analyst  initiates  a  process  by  which  the 
associated  waveforms  are  brought  up  in  the  standard  ARS  display  with  the  data 
time  aligned  on  the  predicted  pP  times  corresponding  to  the  selected  trial  depth. 
After  this,  the  data  can  then  be  further  processed  by  the  analyst  using  the  bandpass 
filter  routine  or  any  of  the  other  signal  processing  algorithms  available  in  the  ARS. 
An  ARS  display  for  this  Lop  Nor  earthquake  is  shown  for  selected  stations  in 
Figure  2-18,  where  it  can  be  seen  that  there  are  multiple  candidate  pP  and  sP 
arrivals  at  the  arrival  times  predicted  for  that  selected  trial  depth.  Moreover,  in  this 
case  the  variation  of  the  pP/P  amplitude  ratios  over  wide  ranges  in  distance  and 
azimuth  permit  the  analyst  to  identify  pP  and  sP  arrivals  with  high  confidence, 
thereby  enabling  this  event  to  be  confidently  screened  out  as  a  natural  earthquake. 
This  implementation  of  the  algorithm  has  been  employed  by  ARS  analysts  in  their 
review  of  the  data  recorded  from  a  large  number  of  explosions  and  shallow 
earthquakes  near  Lop  Nor  as  part  of  the  Advanced  Concept  Demonstration  (ACD), 
which  was  focused  on  improving  the  capability  for  monitoring  the  Chinese  nuclear 
test  site  at  Lop  Nor.  As  part  of  that  demonstration  a  total  of  205  events  were 
thoroughly  analyzed  by  trained  analysts.  The  results  of  this  large  scale  test 
indicates  that  the  new  algorithm  should  be  considered  for  incorporation  into 
routine  seismic  analysis  procedures  such  as  those  employed  at  the  U.S.  NDC  or  the 
IDC. 
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Figure  2-16.  Network  detection  stack  of  pP  for  the  Lop  Nor  region  earthquake 
of  1999/01/27. 
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Figure  2-17.  Network  detection  stacks  of  pP  (top),  sP  (middle)  and  their 
sum  (bottom)  for  the  Lop  Nor  region  earthquake  of  1999/01/27. 
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earthquake  of  1999/01/27,  shown  time  aligned  on  the  predicted  pP  arrival  times 
corresponding  to  the  selected  candidate  depth  (i.e.  18  km)  from  Figure  2-17. 


SECTION  3 

VALIDATION  OF  CANDIDATE  DEPTH  PHASES 


The  network  beamforming  algorithm  described  above  is  proving  to  be  a  powerful 
tool  for  identifying  candidate  depth  phases.  However,  the  presence  of  a  prominent 
peak  on  the  network  detection  stack  does  not,  by  itself,  constitute  a  positive 
identification  which  is  confident  enough  to  be  used  for  purposes  of  event 
screening.  Rather,  it  identifies  possible  depth  phases  which  would  have  to  be 
validated  by  further  review  and  testing.  That  is,  since  the  positive  identification  of 
a  depth  phase  with  a  time  delay  of  greater  than  about  3  seconds  with  respect  to  P 
conclusively  identifies  an  event  as  being  too  deep  to  be  an  explosion,  it  is 
imperative  that  any  such  identification  be  made  with  a  very  high  level  of 
confidence.  As  a  result,  the  development,  implementation  and  testing  of  improved 
analyst  tools  and  quantitative  tests  will  play  a  key  role  in  the  completion  of  any 
new  depth  phase  identification  system.  One  such  tool,  the  Pearce  (1977,1980)  focal 
mechanism  algorithm,  uses  the  observed  relative  amplitudes  of  pP  and/or  sP  depth 
phases  with  respect  to  the  P  wave  amplitude  at  each  station  in  the  network.  The 
solution  space  of  all  double  couple  seismic  sources  is  searched  for  possible 
solutions  that  are  consistent  with  the  relative  amplitudes  that  are  measured  at  each 
station.  If  no  such  solution  can  be  found,  then  it  must  be  concluded  that  either  the 
candidate  depth  phases  were  misidentified  or  that  the  event  was  not  a  double 
couple  (i.e.  it  must  be  considered  as  a  potential  explosion).  Conversely,  if  a  small 
number  of  consistent  solutions  are  found  which  are  compatible  with  the  regional 
tectonics,  then  the  confidence  in  the  identification  of  the  candidate  depth  phases 
would  be  greatly  increased.  Preliminary  results  (Murphy,  et  al,  2000)  suggest  that 
the  Pearce  focal  mechanism  algorithm  may  provide  a  valuable  supplemental 
analyst  tool  for  assessing  candidate  depth  phases  identified  by  the  network  stacking 
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procedure  described  in  the  previous  section.  Another  tool  that  can  be  used  in  the 
characterization  of  network  detection  stack  peaks  is  the  F  detector.  The  F  detector 
(Blandford,  1974)  operates  on  array  station  beam  data  and  provides  an  output  trace 
which  is  approximately  equal  to  the  time  averaged  power  on  the  beam  divided  by 
N  times  the  average  variance  of  the  individual  channels  and  the  beam,  where  N  is 
the  number  of  array  stations.  F  can  be  approximated  by: 


^  f=i/=i  t=i 


(3.1) 


where  u(t)  represents  the  beam  output  trace  at  time  t  from  an  N-element  array  of 
traces  Ui(t),  and  t  =  1,...,M  represents  an  averaging  window  of  M  samples. 
According  to  Blandford  (1974),  when  a  signal  arrives  with  similar  vector-slowness 
to  the  beam,  F  increases  because,  (1)  the  beam  power  in  the  numerator  increases, 
and  (2)  the  denominator  is  reduced  to  the  residual  noise.  If  a  signal  arrives  from  a 
significantly  different  vector-slowness  to  the  beam,  F  is  reduced  as  the 
denominator  increases  faster  than  the  numerator.  The  potential  value  of  the  F 
detector  in  this  context  is  that  it  provides  a  quantitative,  statistical  means  of 
assessing  the  likelihood  that  a  given  arrival  has  the  same  vector  slowness  as  the 
beam.  Therefore,  it  provides  a  potential  mechanism  for  evaluating  whether 
detections  at  a  station  have  characteristics  that  are  consistent  with  what  would  be 
expected  for  pP  and  sP  arrivals  from  the  origin.  This  could  provide  a  powerful 
means  for  eliminating  false  candidate  peaks  associated  with  source  complexity  or 
“mystery  phases”  like  those  identified  in  Section  2.  We  have  implemented  this 
algorithm  in  order  to  evaluate  it’s  capabilities  as  a  tool  to  assist  in  characterizing 
our  depth  phase  stack  peaks.  For  example.  Figure  3-1  shows  the  results  for  an 
earthquake  on  July  9, 1997  located  in  the  Honshu  region  of  Japan  and  recorded  at 
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Figure  3-1.  Best  beam  trace  (top),  F-trace  (middle),  and  probability  trace  (bottom) 
at  station  YKA  for  the  Central  Honshu  earthquake  of  1997/07/09,  mb  =  3.94.  The 
depth  phase  predictions  are  based  on  the  REB  depth  of  71  km. 


the  YKA  array  in  Canada.  In  the  top  of  the  figure  is  the  beam  trace  steered  to  the 
location  of  the  event  and  filtered  between  0.5  and  3.0  Hz.  The  middle  trace  is  the  F 
detector  trace,  calculated  using  a  2  second  moving  window,  and  the  bottom  of  the 
figure  shows  the  F  probability  trace,  which  is  derived  from  the  F  detector  trace  and 
indicates  the  probability  that  a  phase  is  arriving  at  the  same  vector-slowness  that 
the  beam  was  steered  to.  Also  marked  on  the  traces  are  the  predicted  pP  and  sP 
arrivals  for  the  estimated  REB  depth  of  71  km  for  this  event.  It  can  be  seen  from 
the  beam  trace  that  there  is  energy  arriving  near  the  time  of  the  predicted  pP  phase, 
energy  that  also  appears  on  the  F  trace.  The  probability  trace  on  the  bottom  shows 
that  there  is  a  high  probability  that  this  represents  a  phase  arriving  with  the  same 
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azimuth/slowness  that  the  beam  was  steered  to,  giving  us  confidence  that  it  might 
represent  a  true  depth  phase.  The  sP  phase  results  are  not  as  clear,  however,  with 
no  distinct  phase  appearing  on  any  of  the  traces  at  the  expected  arrival  time.  It 
might  be  that  this  phase  is  buried  in  the  coda  of  another,  unpredicted  phase  arriving 
a  few  seconds  earlier  than  the  expected  time  for  sP.  We  have  also  applied  this 
algorithm  to  the  station  YKA  recordings  of  the  shallow  Lop  Nor  earthquake  of 
1999/01/27  discussed  in  Section  2.  The  results  of  this  application  are  shown  in 
Figure  3-2.  Marked  on  the  figure  are  the  predicted  pP  and  sP  times  for  the  REB 
assigned  depth  of  18  km,  in  agreement  with  the  peak  of  our  depth  phase  stack  from 
Section  2  (cf.  Figure  2-16).  Again,  we  see  evidence  in  these  results  that  the 


Figure  3-2.  Beam  trace  (top),  F-trace  (middle),  and  probability  trace  (bottom) 
at  station  YKA  for  the  Lop  Nor  earthquake  of  1999/01/27,  =  3.90.  The 

depth  phase  predictions  are  based  on  the  REB  depth  of  18  km. 
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candidate  depth  phases  are  being  generated  by  this  source  and  as  indicated  by  the 
probability  trace  we  could  say  with  some  confidence  that  they  correspond  to  the 
depth  phases  pP  and  sP.  So  it  appears  that  this  F  detector  algorithm  might  be  useful 
in  the  characterization  of  peaks  resulting  from  our  network  stacking  procedure. 

In  order  to  validate  candidate  depth  phases,  we  have  developed  a  scheme  for 
weighing  our  depth  phase  stacking  results  with  the  results  of  the  F  detector 
processing.  In  this  procedure  we  start  with  the  depth  phase  stacking  results,  an 
example  of  which  is  shown  in  Figure  3-3  for  the  Hokkaido  earthquake  of 
1997/07/01.  It  can  be  seen  fi-om  this  figure  that  there  are  a  number  of  prominent 
peaks  on  the  network  detection  stack,  suggesting  a  number  of  candidate  depth 
phases  which  need  to  be  further  evaluated.  From  this  display  we  then  extract  the 
detections  for  the  array  stations  contributing  to  the  peaks;  these  are  shown  as 
functions  of  candidate  depth  in  Figure  3-4.  Next,  the  F  detector  is  run  on  the 
waveforms  at  each  of  the  arrays  to  give  us  the  detection  probabilities  as  a  function 
of  depth,  as  shown  in  Figure  3-5.  These  probabilities  are  then  used  to  weight  the 
original  detection  traces  and  the  depth  stack  is  re-calculated  from  the  weighted 
traces.  These  results  are  shown  in  Figure  3-6,  where  it  can  be  seen  that  a  number  of 
the  secondary  stack  peaks  from  Figure  3-3  have  been  suppressed,  while  the  peak 
corresponding  to  the  REB  estimated  depth  has  been  retained.  Thus,  these  results 
show  that  we  are  able  to  use  the  F  detector  output  to  provide  a  valuable 
supplemental  analyst  tool  for  assessing  candidate  depth  phases  identified  by  the 
network  detection  stacking  procedure.  Additional  encouraging  results  are  shown  in 
Figure  3-7  for  the  central  Honshu  earthquake  of  1997/07/09,  where  again  the 
weighting  procedure  is  successful  at  suppressing  secondary  peaks  while  retaining 
the  main  peak.  Figure  3-8  shows  additional  results  for  the  Hokkaido  earthquake  of 
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1997/03/13.  These  results  are  even  more  convincing  with  the  only  remaining  peak 
being  the  one  corresponding  to  the  estimated  REB  depth. 


In  summary,  these  preliminary  results  indicate  that  the  combined  use  of  the 
network  stacking  algorithm  with  the  F  detector  has  excellent  potential  for  greatly 
improving  the  confidence  in  depth  phase  identification. 
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Figure  3-3.  Network  detection  stack  of  pP  for  the  Hokkaido  earthquake  of  1997/07/01 
The  dashed  vertical  line  coincides  with  the  REB  depth  estimate  of  63  km. 


30 


1 

1 

1 

1 

1 

n 

ILAR 

rt-i 

1 

1 

1 

t 

inr 

CMAR 

rt*n 

t 

1 

1 

n 

^  YKA 

n  n 

~|  WRA 

1 

1 

1 

i 

FINES 

1 

1 

1 

1 

UL 

n 

ASAR 

n  n 

MORES 

1 

1 

1 

1 

n  n 

n  n  pdar 

T- 

1 

\ 

\ 

1 

GERES 

r- 

1 

1 

1 

1 

TXAR 

1 — P— 1 — 1  ■  1  1  ■< — ^ ^ — 1 — ' — ■ — 1 

f, ,  ,  ,  ,  ,  n  ii  ;  .1 .1  . .  I  ,  .  I  . . I 

0  100  200  300 


DEPTH 

Figure  3-4.  Depth  traces  for  the  array  stations  of  Figure  3-3.  The 
dashed  vertical  line  coincides  with  the  REB  depth  estimate  of  63  km. 
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Figure  3-5.  F  probability  traces  for  the  array  stations  of  Figure  3-3.  The 
dashed  vertical  hne  coincides  with  the  REB  depth  estimate  of  63  km. 
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Figure  3-6.  Weighted  network  detection  stack  of  pP  for  the  Hokkaido  earthquake  of 
1997/07/01.  The  dashed  vertical  line  coincides  with  the  REB  depth  estimate  of  63  km. 
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Figure  3-7.  Comparison  of  original  (left)  and  F  detector  weighted  (right)  network 
detection  stacks  of  pP  for  the  Honshu  earthquake  of  1997/07/09.  The  dashed  vertical  line 
coincides  with  the  REB  depth  estimate  of  70  km. 
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Figure  3-8.  Comparison  of  original  (left)  and  F  detector  weighted  (right)  network 
detection  stacks  of  pP  for  the  Hokkaido  earthquake  of  1997/03/13.  The  dashed  vertical 
line  coincides  with  the  REB  depth  estimate  of  128  km. 
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SECTION  4 

FOCAL  DEPTH  ESTIMATION  BASED  ON  REGIONAL  S  MINUS  P  TIMES 


As  was  noted  previously  above,  only  20%  of  the  free  depth  solutions  in  the 
GSETT-3  database  were  based  on  depth  phase  observations.  It  follows  that,  even 
if  depth  phase  identification  can  be  significantly  improved  by  employing  new 
techniques  such  as  the  network  beamforming  procedures  described  in  Section  2, 
there  will  still  be  many  earthquakes  for  which  focal  depth  can  only  be  estimated 
through  analyses  of  the  observed  P  wave  first  arrival  times.  Therefore,  it  is 
appropriate  to  consider  the  various  factors  affecting  the  accuracy  of  such  estimates. 
One  obvious,  fundamental  limitation  in  depth  resolution  is  associated  with  the 
magnitudes  of  the  differences  in  the  expected  P  wave  first  arrival  times  between 
events  with  different  focal  depths,  relative  to  the  errors  in  the  arrival  time  readings. 
This  issue  is  graphically  illustrated  in  Figure  4-1,  where  the  differences  in  P  wave 
arrival  times  predicted  by  the  lASPEI  91  tables  between  an  event  with  a  focal 
depth  of  h=0  and  events  with  the  same  epicenter  having  focal  depths  of  100  and 
200  km  are  displayed  as  functions  of  epicentral  distance.  Considering  first  the 
event  with  a  depth  of  100  km,  it  is  interesting  to  note  that  for  distances  greater  than 
about  25®,  the  differences  in  arrival  times  can  be  very  closely  approximated  by  a 
linear  relation  of  the  form  Ax  =  Atq  -  txA  where  Ato*=  -9.9  seconds  and  a  « -0.035 
seconds/degree.  Now  consider  the  hypothetical  situation  where  arrival  time  data 
are  available  from  a  100  km  deep  event  from  a  large  network  with  stations  well- 
distributed  with  respect  to  azimuth  over  the  distance  range  25®  <  A  <  90°. 
Assuming  these  data  contain  no  systematic  errors,  a  hypocenter  solution 
determined  under  the  constraint  h  =  0  will  give  the  correct  epicenter  with  an 
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Figure  4-1.  Differences  in  P  wave  arrival  times  as  a  function  of  distance 
predicted  by  the  lASPEI  91  tables  between  an  event  with  a  focal  depth  of  h  =  0 
and  events  at  the  same  epicenter  having  focal  depths  of  100  and  200  km. 


associated  origin  time  which  is  about  12  seconds  earlier  than  the  true  origin  time. 
This  is  the  often  referred  to  trade  off  between  focal  depth  and  origin  time.  The 
travel  time  residuals  with  respect  to  this  constrained  solution  would  be  uniformly 
distributed  over  the  range  extending  from  —1.0  to  +  1.0  seconds,  which  is 
comparable  to  the  expected  error  bounds  on  observed  teleseismie  P  wave  arrival 
time  data.  Thus,  even  with  station-corrected  data,  teleseismie  P  wave  first  arrival 
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times  for  events  with  focal  depths  of  up  to  100  km  or  more  can’t  be  used  to 
confidently  conclude  that  such  events  are  too  deep  to  be  explosions,  even  with  an 
ideal  network  of  stations  located  at  distances  of  greater  than  25°.  Although  the 
detailed  nature  of  this  trade-off  between  focal  depth  and  origin  time  can  be 
complex  in  practical  applications  to  small  events,  an  ability  to  independently 
constrain  the  event  origin  time  can  nonetheless  be  shown  to  lead  to  significantly 
improved  focal  depth  estimates.  More  specifically,  if  S-P  differential  arrival  time 
information  is  available  from  well-calibrated  regional  stations,  it  can  be  used  to 
provide  useful,  independent  estimates  of  the  event  origin  time  for  use  as  constraints 
in  the  hypocenter  location  process.  That  is,  if  to  is  the  event  origin  time  and  tp  and 
ts  are  the  observed  arrival  times  of  P  and  S  waves  that  have  traveled  essentially  the 
same  path  to  a  station  at  distance  A  from  the  event,  it  follows  that: 


tp  =  to  +T,(A,h) 
h  =  h  +T,{A,h) 


(4.1) 


where  h  denotes  the  focal  depth  and  Tp,  Ts  are  the  travel-time  functions  for  P  and  S 
arrivals  as  a  function  of  epicentral  distance  and  focal  depth.  It  follows  that  we  can 
relate  the  origin  time  to  the  observed  S-P  time  interval  as 


tp  -to 


'TsM 


n(^S  “^P  ) 


(4.2) 


or 


h  - 


nOs  ) 


(4.3) 
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Now,  if  the  travel  paths  for  the  observed  P  and  S  arrivals  are  the  same  and  the  ratio 
of  compressional  (a)  to  shear  {P)  wave  velocity  along  these  paths  is  constant, 
independent  of  A  and  h,  then  (3)  can  be  simplified  to 


^0 


-^-1 

J 


n(^S  -^p) 


(4.4) 


which  is  the  relation  originally  proposed  by  Wadati  some  75  years  ago  in  his 
studies  of  deep,  Japanese  earthquakes.  It  follows  that  if  the  average  ratio  a ip 
along  paths  to  a  selected  regional  station  from  events  in  a  specified  source  region 
can  be  determined,  then  the  event  origin  time  can  be  estimated  from  observed  P 
and  S  arrival  times. 

More  formally,  for  a  fixed  event  with  m  observations,  the  hypocentral  uncertainty 
due  to  random  picking  errors  is  given  by  the  3  x  3  variance  matrix  of  the 
hypocenter,  which  in  turn  is  determined  by  the  slowness  vectors  (at  the 
hypocenter)  for  the  m  rays  that  were  observed.  The  trade-off  of  location  with 
origin  time  is  accounted  for  by  removing  the  mean  of  these  slowness  vectors  prior 
to  computing  the  variance  matrix.  If  all  the  data  are  teleseismic  P  waves,  the 
vertical  components  of  all  the  slowness  vectors  are  approximately  the  same  and 
depth  resolution  will  be  poor.  To  increase  depth  information,  it  is  necessary  to  add 
data  for  which  the  vertical  components  of  the  slowness  vectors  are  significantly 
different  from  those  of  the  teleseismic  P  data.  One  way  to  accomplish  this  is  to  add 
a  regional  S-P  based  origin  time  estimate,  which  is  associated  with  a  slowness 
vector  of  zero.  That  is,  if  the  observed  P  and  S  arrival  times  at  some  regional 
station  are  represented  by  equations  4.1  above,  then  since  tp-to  is  approximately  a 
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linear  function  of  ts-tp,  it  follows  that  Ts  can  be  represented  as  a  linear  function  of 

Tp 


T^{A,h)=a  +  bT^(A.h) 


(4.5) 


where  a,b  are  calibration  constants.  Then,  given  observations  of  to,  ts,  we  can  solve 
for  the  origin  time  as 


t 


o 


a  +  bt 


b-1 


(4.6) 


Note  that,  as  a  consequence  of  the  assumed  linear  relation  between  Tp  and  Ts,  this 
origin  time  estimate  is  independent  of  any  estimate  of  the  hypocenter.  We  are 
currently  developing  a  formal  error  model  which  will  rigorously  account  for  the 
effects  of  uncertainties  in  the  calibration  constants  a  and  b  on  the  corresponding 
uncertainty  in  the  resulting  focal  depth  estimate. 

An  illustration  of  the  power  of  calibrated  regional  S-P  data  in  constraining  event 
focal  depth  can  be  seen  in  Figure  4-2,  where  we  show  a  comparison  of  ts  -  tp 
variation  as  a  function  of  focal  depth  predicted  at  a  fixed  near-regional  distance  of 
5°  by  the  IASPEI91  travel  time  tables  with  the  corresponding  predicted  variation  in 
P  wave  first  arrival  times  at  a  nominal  distance  of  45®.  For  20  <  h  <  100  km,  ts-  tp 
is  essentially  independent  of  depth  in  the  regional  distance  range  2°  <  A  <  10°, 
while  tp  at  A  —  45°  decreases  by  more  than  1  sec  for  every  10  km  increase  in  h. 
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Figure  4-2.  Comparison  of  ts-tp  variation  as  a  function  of  focal  depth  predicted  at 
a  fixed  near-regional  distance  of  5®  by  the  IASPEI91  travel  time  tables  with  the 
corresponding  predicted  variation  in  P  wave  first  arrival  times  at  a  nominal 
distance  of  45°. 

Therefore,  such  regional  data,  if  properly  calibrated,  provide  very  robust 
constraints  on  origin  time  and  the  associated  depth  estimate. 

An  application  of  this  constrained  origin  time  procedure  is  illustrated  in  Figure  4-3, 
where  we  show  a  comparison  of  differences  in  depths  estimated  from  teleseismic  P 
arrival  times  with  respect  to  the  “true”  depth  phase  constrained  values  for  47 
Honshu  earthquakes  in  the  depth  range  of  50  to  150  km,  obtained  with  and  without 
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origin  time  constraints.  These  results  indicate  that  earthquakes  with  depths  greater 
than  about  30  km  in  this  region  can  be  shown  to  be  too  deep  to  be  explosions  at  a 
very  high  level  of  confidence  if  S  -  P  constrained  origin  times  are  available. 


o  Depth  From  Unconstrained  Teleseismic  P  :  Ah  =  1 1 .7  km,  a  =  32.3  km 
•  Depth  From  Origin  Time  Constrained  Teleseismic  P  :  Ah  =  3.7  km,  a  =  8.0  km 


Figure  4-3.  Comparison  of  differences  in  estimated  depth  with  respect  to  the 
“true”  depth  phase  constrained  values  for  47  Honshu  earthquakes  in  the  depth 
range  of  50  to  150  km  obtained  with  and  without  origin  time  constraints. 


We  have  formulated  a  focal  depth  hypothesis  test  by  comparing  this  S-P-based 
origin  time  estimate,  fo(S-P),  with  the  corresponding  teleseismic-based  origin  time, 
^  (ft = 0  km),  determined  with  the  depth  constrained  to  zero.  Substantial  differences 

in  these  two  estimates  provide  evidence  of  significant  focal  depth.  In  particular, 
since  from  the  IASPEI-91  tables  ^(ft  =  10  km)-  ^(ft  =  0  km )«  1.5  seconds  and  the 
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measured  variability  of  (s  -  p)-  ^  (/i  =  0  km )  for  well-located  events  recorded  at 

calibrated  regional  stations  is  on  the  order  of  ±  3.5  seconds,  we  have  assumed  for 
preliminary  analysis  purposes  that  differences  of  (s  -  P)-  ^  (A  =  0  km  )  greater 

than  about  5  seconds  provide  high  confidence  that  the  event  depth  is  greater  than 
10  km  and,  consequently,  that  the  event  can  be  confidently  screened  out  as  a 
natural  earthquake  source. 

We  are  currently  evaluating  this  hypothesis  test  using  regional  data  recorded  from 
events  in  the  source  regions  of  Figure  2-11  for  which  the  event  depths  have  been 
well  constrained  by  validated  depth  phase  observations.  For  this  application,  P  and 
S  arrival  times  were  determined  at  selected  regional  stations  for  well  located  REB 
events  and  the  resulting  t?  -  to  and  ts  - time  intervals  were  used  to  evaluate  the 
calibration  constant  in  equation  4-4,  where  to  denotes  the  origin  time  from  the  well 
constrained  hypocenter  solution.  Origin  times  for  events  recorded  at  each 
calibrated  station  were  then  estimated  from  observed  ts  -  h  intervals  and  compared 
with  corresponding  constrained  teleseismic  origin  time  estimates,  ^  (/i=okm),  to 

test  for  significant  focal  depth.  This  process  is  illustrated  for  a  sample  of  Honshu 
events  recorded  at  station  MJAR  in  Figure  4-4,  where  the  determination  of  the 
calibration  constant  (left)  and  its  application  in  a  focal  depth  hypothesis  test  (right) 
are  shown.  It  can  be  seen  that,  in  this  case,  almost  all  the  events  with  validated 
focal  depths  of  greater  than  50  km  can  be  screened  out  as  being  too  deep  to  be 
explosions  with  high  confidence.  The  dashed  line  in  the  right-hand  panel  of  this 
figure  denotes  the  value  of  ^(/j)  _  (ft = Okm)  predicted  as  a  function  of  focal  depth 

by  the  teleseismic  IASPEI-91  travel-time  values,  and  it  can  be  seen  that  the 
observed  data  show  some  systematic  deviations  from  this  relation  that  may  indicate 


41 
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h,  km 


Figure  4-4.  Application  of  the  S-P  origin  time  hypothesis  test  to  data  recorded  from 
selected  Honshu  events  at  regional  station  MJAR  (A  =  1.9°). 


variations  in  the  regional  station  calibration  constant  as  a  function  of  epicentral 
distance  and  focal  depth  that  will  need  to  be  further  evaluated.  In  any  case,  this 
hypothesis  test  clearly  shows  promise  for  applications  to  event  screening  based  on 
focal  depth.  The  corresponding  calibration  results  for  the  Hokkaido  region  of  Japan 
are  shown  in  Figure  4-5  where  again,  almost  all  events  with  validated  depths  of 
greater  than  50  km  are  screened  out  as  being  to  deep  to  be  explosions  with  high 
confidence.  In  this  case,  the  Hokkaido  points  in  the  right  hand  plot  show  larger 
scatter  relative  to  the  lASPEI  prediction.  This  could  be  due  to  variation  in  the 
paths  sampled  or  errors  in  the  S  or  P  arrival  time  picks,  areas  currently  under 
investigation.  Finally,  the  calibration  results  for  station  NIL  for  events  in  the  Hindu 
Kush  region  of  Asia  are  shown  in  Figure  4-6.  Once  again,  all  events  are  screened 
out  as  being  to  deep  to  be  explosions,  and  here  the  data  in  the  right-hand  plot  more 
closely  follow  the  Uend  predicted  by  the  simple,  lASPEI  earth  model. 
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Figure  4-5.  Application  of  the  S-P  origin  time  hypotl^sis  test  to  data  recorded  from 
selected  Hokkaido  events  at  regional  station  MJAR  (A  ~  7.6°). 
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Figure  4-6.  Application  of  the  S-P  origin  time  hypothesis  test  to  data  recorded  from 
selected  Hindu  Kush  events  at  regional  station  NIL  (A  »  4.7°). 
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The  stability  of  the  S-P  calibration  procedure  can  be  evaluated  by  examining 
events  for  which  we  have  obtained  absolute,  independent  knowledge  of  the  event 
depths  and  origin  times.  In  this  regard,  we  have  collected  P  and  S  arrival  time 
information  from  a  large  sample  of  PNE  explosions  recorded  at  station  BRV, 
where  the  event  depths  and  origin  times  are  accurately  known.  Figure  4-7  shows 
the  calibration  curve  (left)  and  associated  errors  in  S-P  origin  time  estimates 
(right),  where  these  errors  are  calculated  with  respect  to  the  known  explosion 
origin  times.  Here  we  can  see  that  the  explosion  data  appear  to  show  larger  scatter 
than  that  seen  in  the  earthquake  data.  This  could  be  due  to  a  greater  variation  in  the 
travel  paths  for  the  explosions  relative  to  the  earthquakes  for  these  widely 
dispersed  PNE  events,  as  well  as  the  generally  weak  S  wave  onsets  observed  from 
explosion  sources. 


(ts  -  tp),  sec  A,  degrees 


Figure  4-7.  Calibration  curve  (left)  and  associated  errors  in  S-P  origin  time  estimates 
(right)  obtained  for  Soviet  PNE  events  recorded  at  Borovoye. 
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Figure  4-8  shows  the  errors  in  the  S-P  origin  time  estimates  as  a  function  of  event 
location  around  Borovoye.  It  can  be  seen  here  that  the  travel  paths  sampled  by 
these  Borovoye  recordings  cover  a  broad  geographic  area.  However,  this  figure 
also  shows  that  there  does  not  appear  to  be  any  correlation  between  the  size  and 
sign  of  the  origin  time  residuals  and  the  event  locations. 

-5  sec  +5  sec 

OOOoo»*++H — I - 


Figure  4-8.  S-P  origin  time  residuals  as  a  function  of  event  location  for  PNE’s 
recorded  at  station  BRV. 
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Clearly,  this  method  of  independently  constraining  the  event  origin  time  shows 
promise,  in  that  it  enables  us  to  screen  out  events  as  having  been  too  deep  to  be 
explosions  in  the  absence  of  any  depth  phase  information.  It  might  also  be  used  as 
a  supplemental  check  on  the  validity  of  depth  phase  picks  when  the  move-out  of 
the  picks  is  small  or  there  is  large  uncertainty  in  our  confidence  of  the  depth  phase 
picks.  We  are  continuing  to  collect  additional  regional  data  from  well  located 
events  in  an  attempt  to  evaluate  the  reliability  of  this  origin  time  constraint  and  to 
determine  its  dependence  on  source  location  and  depth.  We  are  also  investigating 
how  best  to  incorporate  this  S-P  calibration  information  into  a  more  formal 
computation  of  focal  depth  uncertainty. 
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SECTIONS 

SUMMARY  AND  CONCLUSIONS 


This  report  has  provided  a  summary  of  the  results  of  a  research  program  that  has 
been  directed  toward  the  development  of  improved  procedures  for  estimating  focal 
depths  of  seismic  sources  for  use  in  the  monitoring  of  underground  nuclear 
explosions.  In  particular,  new  analyst  tools  have  been  implemented  and  tested  in 
an  attempt  to  increase  the  detection  and  identification  of  the  depth  phases  pP  and 
sP,  and  a  method  of  obtaining  independent  estimates  of  event  origin  time  through 
the  use  of  differential  S-P  arrival  times  has  been  investigated. 

The  detection  and  identification  of  the  secondary  depth  phases  pP  and  sP  were 
addressed  in  Section  2,  where  a  number  of  examples  were  presented  to  illustrate 
the  effects  of  complicating  factors  such  as  the  relative  frequency  content  of  pP  and 
sP  with  respect  to  P,  earthquake  focal  mechanism  and  the  occurrence  of  anomalous 
P  phases  which  are  not  predicted  by  the  global  travel  time  tables.  A  fully 
automatic  network  beamforming  algorithm  was  then  described  for  stacking  raw 
IMS  detection  data  to  more  consistently  identify  candidate  depth  phases  for 
subsequent  review  by  the  IDC  analyst.  This  new  procedure  has  been  tested  and 
evaluated  using  data  recorded  from  selected  events  in  the  Pamir-Hindu  Kush, 
Hokkaido  and  central  Honshu  seismic  zones  as  well  as  on  events  located  near  the 
Chinese  explosion  test  site  at  Lop  Nor.  This  algorithm  has  been  integrated  into  the 
Analyst  Review  Station  which  is  employed  at  both  the  U.S.  NDC  and  the  IDC  and 
can  now  be  tested  using  the  full-processing  environment  of  the  representative 
nuclear  monitoring  systems. 
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In  Section  3  we  investigated  the  usefulness  of  the  F  detector  and  other  tools  for 
validating  candidate  depth  phases  identified  by  our  network  beamforming 
algorithm.  That  is,  the  presence  of  a  prominent  peak  on  a  network  detection  stack 
does  not,  by  itself,  constitute  a  positive  identification  which  is  confident  enough  to 
be  used  for  purposes  of  event  identification.  Rather,  it  identifies  possible  depth 
phases  which  have  to  be  validated  by  further  review  and  testing.  We  showed  that 
the  development,  implementation  and  testing  of  improved  analyst  tools  and 
quantitative  tests  will  play  a  key  role  in  the  development  of  any  new  depth  phase 
identification  system. 

A  brief  description  of  focal  depth  estimation  determined  from  P  wave  first  arrival 
times  alone  was  presented  in  Section  4,  where  we  describe  a  method  of 
determining  independent  event  origin  times  using  regional  S  and  P  arrival  time 
information.  Specifically,  it  was  shown  that  if  S-P  differential  arrival  time 
information  is  available  from  well-calibrated  regional  stations,  it  can  be  used  to 
provide  useful,  independent  estimates  of  the  event  origin  time  for  use  as  constraints 
in  the  hypocenter  location  process.  This  method  was  applied  to  a  series  of  well- 
located  events  from  the  Hindu  Kush  region  and  the  Hokkaido  and  Honshu  regions 
of  Japan  with  promising  results. 

The  research  summarized  above  supports  the  following  conclusions  regarding 
seismic  determination  of  focal  depth: 

(1)  Analyst  detection  and  identification  of  the  depth  phases  pP  and  sP  is 

complicated  by  the  effects  of  earthquake  focal  mechanism  and  a  variety  of 
other  factors.  As  a  result,  depth  phase  data  were  reported  for  only  about  5%  of 
the  REB  events  analyzed  during  the  GSETT-3  experiment. 
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(2)  A  new,  automatic  processing  algorithm  based  on  network  beamforming  of  raw 
IMS  detection  data  has  been  developed  which  has  been  shown  to  identify 
significant  numbers  of  candidate  depth  phases  which  were  not  detected  by  the 
PIDC  analysts.  Application  of  this  algorithm  to  data  from  a  large  sample  of 
earthquakes  from  the  Lop  Nor,  Pamir-Hindu  Kush,  Hokkaido  and  central 
Honshu  seismic  zones  resulted  in  the  identification  of  prominent  candidate 
depth  phase  peaks  near  the  published  REB  depth  for  over  70%  of  the  events 
analyzed,  including  some  with  magnitudes  as  low  as  mb=3.75  and  depths  as 
shallow  as  18  km.  This  algorithm  has  been  integrated  into  the  Analyst  Review 
Station  (ARS)  software  which  is  employed  at  both  the  U.S.  NDC  and  the  IDC, 
and  can  now  be  evaluated  using  the  full-processing  environment  of  the 
representative  nuclear  monitoring  systems  as  well  as  historical  data  archives. 

(3)  Preliminary  applications  of  the  F  detector  algorithm  to  the  confirmation  of 
depth  phase  identification  have  been  encouraging  and  suggest  that  it  could  be 
used  in  conjunction  with  the  network  beamforming  algorithm  to  provide 
analysts  with  a  much  improved  capability  to  confidently  identify  more  pP  and 
sP  arrivals. 

(4)  Analyses  of  focal  depth  estimates  determined  using  P  wave  first  arrival  times 
alone  indicate  that  there  are  a  number  of  sources  of  systematic  error  which  are 
not  accounted  for  by  existing  hypocentral  inversion  codes  and  which  can  lead 
to  significant  bias  in  estimated  focal  depths  and  their  associated  uncertainty 
bounds.  This  shortcoming  seriously  limits  the  current  applicability  of 
estimated  focal  depth  to  nuclear  monitoring.  However,  we  have  shown  that  if 
S-P  differential  arrival  time  information  is  available  from  well-calibrated 
regional  stations,  it  can  be  used  to  provide  useful,  independent  estimates  of  the 
event  origin  time  for  use  as  constraints  in  the  hypocenter  location  process. 

That  is,  this  method  of  independently  constraining  the  event  origin  time  could 
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be  used  to  confidently  screen  out  events  as  having  been  too  deep  to  be 
explosions  in  the  absence  of  any  depth  phase  information. 
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